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ABSTRACT 

We formulate the Lagrangian perturbation theory to solve the non-linear dynamics of 
self-gravitating fluid within the framework of the post-Newtonian approximation in 
general relativity, using the (3+1) formalism. Our formulation coincides with Newto- 
nian Lagrangian perturbation theory developed by Buchert for the scale much smaller 
than the horizon scale, and with the gauge invariant linearized theory in the longitu- 
dinal gauge conditions for the linear regime. These are achieved by using the gauge 
04 ' invariant quantities at the initial time when the linearized theory is valid enough. The 

post-Newtonian corrections in the solution of the trajectory field of fluid elements are 
calculated in the explicit forms. Thus our formulation allows us to investigate the evo- 
, lution of the large-scale fluctuations involving relativistic corrections from the early 

regime such as the decoupling time of matter and radiation until today. As a result, 
we are able to show that naive Newtonian cosmology to the structure formation will 
' be a good approximation even for the perturbations with scales not only inside but 

also beyond the present horizon scale in the longitudinal coordinates. Although the 
■ post-Newtonian corrections are small, it is shown that they have a growing trans- 

verse mode which is not present in Newtonian theory as well as in the gauge invariant 
Q ■ linearized theory. Such post-Newtonian order effects might produce characteristic ap- 

pearances of the large-scale structure formation, for example, through the observation 
of anisotropies of the cosmic microwave background radiation (CMB). Furthermore 
since our approach has a straightforward Newtonian limit, it will be also convenient 
for numerical implementation based on the presently available Newtonian simulation. 
Our results easily allow us to perform a simple order estimation of each term in the 
jJJ ■ solution, which indicates that the post-Newtonian corrections may not be neglected in 

the early evolution of the density fluctuation compared with Newtonian perturbation 
solutions. 

Key words: gravitation-instability-cosmology: theory-large-scale structure of Uni- 
verse 



1 INTRODUCTION 



The observation of isot ropy of the cosmi c microwave background (CMB) indicates that the universe is remarkably isotropic 
over the horizon scale ( Smoot et al.1992 ). Thus it is natural to describe the large scale spatial geometry of the universe by 
a homogeneous and isotropic metric, namely, the Friedmann-Robertson-Walker (FRW) model. However, the real universe is 
neither isotropic nor homogeneous on local scales and has a hierarchical structure such as galaxies, clusters of galaxies, super- 
clusters of galaxies and so on. In the standard big bang scenario, it is commonly believed that such large-scale inhomogeneous 
structures of the universe were formed as a result of growt h of d ensity fluctuations whose amplitudes had been very small 
in the early regime in FRW background (see, e.g., Peebles (1980)). Furthermore, such inhomogeneous structures are usually 
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considered to have been formed only due to gravitational instability after matter (baryons) and radiation decoupled. There- 
fore, the gravitational instability of self-gravitating fluid has been a main subject of the research in cosmology in connection 
to the formation of structures. 

The investigation of the nonlinear large-scale structure formation of the universe have been usually studied within the 
Newtonian theory. This is based on the following fact: as far as the scales of the considered density fluctuations are much 
smaller than the horizon scale, the Newtonian approximation will be accurate enough to describe their evolutions up to 



non- linear regime where the density contrast becomes larger than unity (Peeblesl980). Actually, many large-scale simulations 



based on the Newtonian approximation have been performed to compare the numerical results with the observed distributions 
of galaxies. Analytical approaches to the evolution of density fluctuation based on the Newtonian Lagrangian picture have also 
been developed by Buchert ( 198S| : 



1992 



1994) and the other authors (Bouchet et al.1995; Bouchet et al.1995 



Catelanl995 



Moutarde et al.1991) where it is expected that the approach gives a good approximation up to a certain stage of non- linear 



regime (Buchert, Mellot fc Weifil994; Mellot, Buchert fc WeiB1995; Wei6,Gottlober fc Buchert 1996). Buchert (1992; 199c) 



showed that the general first-order solutions for an Einstein-de Sitter background include well-known Zel'dovich approximation 



( Zel'dovichl970 ) as a special case. 



On the other hand, the region of observations of the large-scale structure is steadily increasing. For example,the very 
ambitious observational program Sloan Digital Sky Survey (SDSS; e.g., Gunn & Weinbergl99E ) will complete a spectroscopic 
survey of ~ 10 6 galaxies brighter than r' ~ 18 mag over it steradians, which covers the region of several hundred megaparsecs. 
As mentioned above, Newtonian approximation is valid only for system with scales much smaller than the horizon scale, so 
it is not clear at all if the application of the Newtonian theory is appropriate for such a wide region of spacetime. In fact the 
fluctuations relevant for the large-scale structures have not been always much smaller than horizon scales in the past. For 
example, the horizon scale at the decoupling time is of the order of rfd cc (l + Zdcc) ~ 80/t -1 Mpc in the present physical length. 
This suggests that we have to employ the relativistic description for the evolution of fluctuations larger than or equivalent to 
such scales. 

Thus to understand the evolution of the large-scale structure of the universe closely, it is important and necessary to have 
some formalism to evaluate the effect of general relativistic corrections to the Newtonian dynamics on such a region. Conversely, 
until we obtain such a formalism, we are not confident about the validity of the Newtonian cosmology in studying evolutions of 
perturbations. The purpose of the present paper is to provide such a formalism based on the Lagrangian perturbation theory 



in the post-Newtonian (PN) framework. In the PN approximation to cosmology (Futamasel996; Shibata & Asadal995), it 
has been shown that the metric is given in the following form at the Newtonian limit 



ds = - 1 



(cdt) 



2 . 2 

+ a 



® (i 



and the Newtonian-like gravitational potential <j> is generated by the density fluctuation field 8 via the cosmological Poisson 
equation = AirGa 2 ptS, where p^ denotes the background density. This includes the non-linear situations for the density 
contrast field where is much larger than unity. The above metric is also applied in deriving the cosmological lens equation 



for a realistic inhomogeneous universe (Futamascl995; Pyne & Birkinshawl996). For recent years, the studies of the gravita- 



tional lensing in the cosmological circumstance have been widely developed using the equation. For example, Seljak (199f) 



investigated the secondary effect of gravitational lensing on CMB anisotropics induced by the non-linear large-scale struc- 
tures using the above Newtonian metric (longitudinal gauge). The future generation of detectors of CMB anisotropies with 
great accuracy requires theoretical efforts to obtain more precise predictions for the anisotropies. In such a case we need 
to consider the gravitational lensing of photons traveling across the present horizon scale (~ 3000ft - 1 Mpc). This indicates 
that, if using the above metric, we first have to investigate the validity whether the cosmological Newtonian metric is able 
to describe accurately the dynamics of fluctuations within such a scale. The following considerations may explain our results 
of the present paper intuitively. The cosmological Poisson equation suggests that the metric perturbation quantity (f>/c 2 is 
of the order of (al / '(ct)) 2 5(n for the density contrast 8m with the characteristic scale I in the present physical length. Then 
we expect that the first PN metric perturbations which determine the first PN force of equations of motion are of the order 
of (0/c 2 ) 2 , namely, (al/ (ct)) 4 8f l} . This means that the quantity (al/(ct)) 8m is an essential expansion parameter in the PN 
approach to cosmology. As mentioned later in detail, the fact that the quantity (I /(do)) 2 8 (to) is much smaller than unity 
for the present horizon-scale fluctuation with the usual power spectra explains why the PN approximation is sufficiently able 
to describe the evolution of the fluctuations not only inside but also beyond the present horizon scale do- Furthermore, the 
formalism based on the PN approximation allows us, for example, a more accurate treatment of propagation of photon on 
the geometry produced by matter inhomogeneities, as required in the study of gravitational lensing by cosmic structures (see, 
e.g., Schneider, Ehlers, & Falcol992) and other applications. 



Naturally, we also re quire our formalism to agree with the gauge invariant linearized theory developed by Bardeen ( 188j ) 
and Kodama & Sasaki ( 1984 ) in the linear regime. There have been some st udies to inv e stigate the nonlinear large-scale 



structure formation based on the relativistic Lagrangian perturbation theory (Kasail993 



Kasail995 



Matarrese, Pantano 



Sac z!994; Matarresel99(: ). However, because of the gauge condition adopted in these studies, namely, the synchronous 
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comoving coordinates, it is not easy to have contact with the Newtonian Lagrangian approach fully developed and used for 
the numerical simulation. In particular, Matarrese & Terranova (1996) called their formalism post-Newtonian treatment in 
spite of starting with the synchronous comoving coordinates. As is well known, a peculiar velocity of any fluid element does 
not appear explicitly in such a coordinate, which makes the PN treatment difficult. Thus we will not restrict ourselves to such 
a coordinates, but we shall use (3+1) formalism of general relativity to admit various choices of the coordinates. In particular 
we shall study the Lagrangian perturbation theory in the coordinates where the comparison with the Newtonian case is most 
easily made. In this way our formalism would be also convenient to use for the numerical simulation straightforwardly based 
on Newtonian simulations. 

There have been another relativistic Lagrangian approaches. Bertschinger & Hamilton (1994) and Bertschinger & Jain 
( 1994 ) have focused on the Lagrangian evolution equations of the electric and magnetic parts of the Weyl tensor for cold dust 
in the Newtonian limit. They derived the equations using not Einstein equations but the relativistic kinematical prescription 
for fluid flow developed by Ehlers (1961) and Ellis (1971). They suggested that the magnetic part does not vanish generally 
in the Newtonian limit, implying that the Lagrangian evolution of the fluid is not purely local. This is a so-called "non-local" 
problem also discussed by the other authors (Hui & Bertschingerl996; Kofman & Pogosyanl995). Here, we formulate another 
relativistic Lagrangian perturbation theory for fluid motion from different point of view. 

Recently, Matarrese, Mollerach, & Bruni (1998) derived the second order solutions of metric perturbations in the Poisson 
gauge that reduces to the longitudinal gauge in the case of scalar first order perturbations. They obtained the solutions by 
transforming the known solutions in the synchronous comoving coordinates to them in the Poisson gauge, using the second 
order gauge transformation. We here start with the (3+1) formalism which includes the Poisson gauge from the beginning . 
Therefore, our PN treatment may be more transparent physically, and can give a physical understanding of the validity of 
the Newtonian treatment. Actually, the time dependence of our second order solutions agrees with their results as described 
below. 

This paper is organized as follows. In Section 2 we shall write down the cosmological PN approximation in (3+1) formalism. 
In Section 3 we introduce the Lagrangian perturbation theory to iteratively solve the PN equations derived in Section 2. In 
Section 4 we summarize our results and discuss the effect of PN corrections and evaluate explicitly its importance. Summary 
and conclusions are given in Section 5. Throughout this paper, Greek and Latin indices take 0, 1, 2, 3 and 1, 2, 3, respectively. 



2 THE POST-NEWTONIAN APPROXIMATION IN (3+1) FORMALISM 



2.1 The (3+1) formalism in cosmological situation 

We wish to develop the relativistic Lagrangian perturbation theory formulated in the coordinates system which has a straight- 
forward Newtonian limit. For this purpose, it is convenient to use th e (3+1) formali sm. The post-Newtonian (PN) approxima- 
tion in the (3 +1) f ormalism has been formulated by one of authors ( Futamasel992) and further developed b y Asada, Shibata 



Shibata fc Asadal995). We first review 



& Futamase ( 1996 ). Its cosmological generalization is straightforward ( Futamasel996 ; 
the (3+1) formalism in the cosmological situation. 

Here we shall present the basic equations using the (3+1) formalism. Let us first assume that there exists a congruence 
of timelike worldlines from which the spacetime looks isotropic. We shall call the family of the worldlines basic observers, who 
see no dipole component of anisotropies of the cosmic microwave background (CMB). We can regard that any one of these 
observers at the spacetime point x moves with 4- velocity n M (a;) without loss of generality. The tangent vector n M is normalized 
as n M n M = —1. These observers are used to foliate the spacetime by their simultaneous surfaces: t = const. We shall assume 
that the three surfaces form hypersurfaces with the unit normal vector field n M . 

In the (3+1) formalism, the metric is split as 



n M = (-a,0), 

1 (3 V 



n 



a 



a 



(2.1) 



where 7^ is the projection tensor and a, /3 l , and 7^ are the lapse function, shift vector and metric on the spatial hypersurface, 
respectively. Then the line element may be written as 



ds 2 = — (n M dx M ) 2 + 'y l _ lv dx^dx v 

= — (q 2 — (3if3 l )(cdt) 2 + 2f3i(cdt)dx 1 + ■y i jdx 1 dx : ' 



(2.2) 



We shall keep the power of c explicitly henceforth because we make the 1/c expansion in various quantities in the post- 
Newtonian approximation. 
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We may define the extrinsic curvature of the hypersurface by 

Kij = — ru;j (2.3) 

where ; denotes the covariant derivative with respect to the metric g^ v . The basic variables in the (3+1) formalism are then 
the spatial metric and the extrinsic curvature. The lapse function and shift vectors are freely specified as the coordinate 
conditions. The (3+1) formalism naturally decomposes the Einstein equations, 

_ 8tyG _ 

into the constraint equations and the evolution equations. A denotes the cosmological constant. The former set of equations 
constitutes the so-called Hamiltonian and momentum constraints. These become 

R- K, 3 K 13 +K 2 = ±^£ + 2A, (2.4) 

DtlCj - DjK = —Jj, (2.5) 

where K — y 1 - 1 Kij, Di and R are the covariant derivative and the scalar curvature with respect to yij, respectively. E and Jj 
are defined as 

J* = -T tu ,ra'Vi. (2.6) 
The evolution equations for the spatial metric and the extrinsic curvature become 

~ JW« = -2aKij + Aft + Aft, (2.7) 
1 d 

cdi Kii = + KKij ~ 2KaK i) ~ D i D i a ~ aA 7v 

+(DjP m )K mi + {D % f3 m )K mj + /3 m D m Kij - ^fa (Sij + \lij{E ~ S'i)j , (2.8) 
where Rij is the Ricci tensor of the hypersurface, and 

Si^WW (2.9) 

For the sake of convenience we shall write down the evolution equations for the trace of the spatial metric and the extrinsic 
curvature which will be sometimes used in the following: 

~T = 2r / (-aK + D i l ), (2.10) 
Id „ ,„ . . „,„ „ 4nG 



cdt 



K = a{R + K*) - D l D r a + ft DjK + -^a(S\ - 3E) - 3aA, (2.11) 



c 



where 7 is the determinant of 7^ . 

For the cosmological application we introduce the scale factor which expresses the expansion of the universe. Following 
Ehlers (1961), we define the spatial distance between two nearby timelike worldlines with tangent vector as 

«=(7m^V) V2 , (2.12) 

where rf is the vector connecting these two worldlines. Then the change of the distance along the worldline may be calculated 
as 

4-(S£) = n M (^) ;M = -[ \k + <TycV]« (2.13) 
dr 3 

where e 1 = rf/5l, with rf = 7 I /J '7 A \ is the unit spatial vector and ct^ is the trace-free part of the extrinsic curvature. It is 
natural to interpret the above equation (2.13) as expressing the cosmic expansion, that is, the first and the second terms on 
the right-hand side represent the isotropic and anisotropic expansion, respectively. Motivated by this interpretation, we shall 
introduce the scale factor as 

K = ---, (2.14) 
c a 

where 

da da 

a{t) =m x =di' 

where a — a(t) is the scale factor as a function of the cosmic time t. This gives a slice condition for a corresponding to the 
maximal slicing in an asymptotically flat case. For more general slice conditions, we shall take the following form for K: 

K=---+6(x). (2.15) 
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It turns out that it is convenient to use the following variables instead of the original variables: 



R = a R. 



We also define cry as 



(2.16) 



(2.17) 



We should note that in our notation, indices of cry are raised and lowered by 7y, so that the relations a j = a l j and a 13 = c?a li 
hold. Using these variables, the constraint equations may be rewritten as 



1 cV 
c 2 a 2 

DiO z i — 



6a 2 
2 
3 



R- -a 

6 



'J 



3c a Q 3<~ 4 



9 



3c J 



(2.18) 
(2.19) 



where A and R are the covariant derivative and the scalar curvature with respect to 7y, respectively. The evolution equations 
take the following forms: 



13- 2 d_ 



2 1 / - 

2acry - -a6»7y + (A/3j + A/% 



:9i 



a \R l 



1 1 9 . 6, 

--^-7 = -(a 

c 7 at c 

3 a 3 , 
-o - + -T (a - 

er a c 2 



- /-') - (ad, 

2ac> li a i J + (Di/3 l )aij + {bjl3 l )a u + I3 l (bia l:j 
2 (-Q0 + A/3 1 ) , 



c 

8ttG 



-(a + 2) 



a , 



+ ■ 



D- 

a 



180 
~c~di 



2d £ 
c a 3 



1 Q ; 

3 7*7-3 1 



^Difl - -^a(S\ + E) + aA, 



(2.20) 

(2.21) 
(2.22) 

(2.23) 



where 7 denotes the determinant of 7y . 

We shall now introduce the spatial background geometry motivated by the following consideration. It seems that the 
metric perturbations in the present universe remain small almost everywhere even when the density contrast is much larger 
than unity. Thus it is natural to assume that the metric structure of the universe may be described by a small perturbation 
from the Friedmann-Robertson- Walker (FRW) spacetime in an appropriate coordinate system. Then, we shall further specify 
the spatial metric as (Bardeenl882; |Kodama fc Sasakil98^ ) 

(2.24) 



7y = (1 - 2^)7y B) + hij, 
where 7 I - S ' is the spatial metric of the FRW geometry: 
1 



0, 



-(B) 

7y 



+ f 



(2.25) 



where /C represents the curvature parameter of FRW models and r 2 = x\ + x\ + x\. We treat the perturbations as being small 
quantities, so we are able to regard spatial tensorial quantities in our equations as tensors with respect to the background 
spatial metric 

Our metric is comp letely general. To reduce the gauge freedom we impose the following transversality constraints on hij 
(|5hibata fc Asadal995|): 



k "»J 



0, 



(2.26) 



denotes covariant derivative with respect to 7L , whose inverse is 7 



-JB)ij 



By means of the condition (2.26), we 



7 

where D 

can erase the scalar part and the vector part in fey, and it guarantees that hij only contains the tensor mode in the PN order 
which represents the freedom of the gravitational wave. Thus we need to take into account only the linear term in hij because 
we are interested in the first PN approximation in the present paper. If we erase the scalar part of the shift vector by using 
the residual gauge freedom, the coordinate condition corresponds to the so-called "longitudinal gauge" (with respect to scalar 
modes) in the linear theory or "Poisson gauge" discussed by Bertschinger (1995) and Ma & Bertschinger (1995). However, we 
leave the freedom here for the present and adopt only the condition (2.26) for the sake of generality. Then the inverse matrix 
of 7y become 

1 



~' l J 

7 = 



1 - 2^ 



{B)ij 



h'- 



(2.27) 



where h 1 j = ; y < - B ^ k hkj- Using the above variables, we can rewrite Rij 



Rij 



2£7lf + - 



2i> 
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+- 



+ 3/C/ii 



i£» (s) £) (s)fe /i 



(2.28) 



(l-2^) 2 

The above discussion is adequate even if the universe has the nonlinear structures. 

As one of the present authors has pointed out ( Futamasel989; Futamasel996 ), the problem of the back reaction of the 
local inhomogeneities on the global expansion arises when the universe has nonlinear structures (8 ^> 1) on small scales. 
However, we are here interested in the evolutions of the fluctuations from sufficiently early time such as the decoupling time 
of matter and radiation. Then the linear theory will be valid in the early evolution. Thus we assume that we can define the 
background to obey Friedmann law until the universe has quasi-nonlinear structures (5 ;> 1). In other words, we assume that 
even if the back reaction for the expansion of background exists, it is negligibly small so that we can perform the perturbative 
approach superimposed on FRW cosmologies. Accordingly, from the lowest order of Eq.(2.18), we first set 

J &nGp b {t) K.c 2 , c 2 A ( .,.„ )1 



Kc 2 

„2 



3 ' 



where H(t) and pb(t) are the Hubble parameter and the homogeneous density of the background FRW universe, respectively. 
For simplicity, we may regard pb as an averaged value over the volume as large as the horizon scale (c£) 3 : (p)f c4 )3 = Pb(t)- The 
continuity equation shows that the averaged mass density is conserved, namely, pbO 3 is constant for the dust model. 

In the following, we only consider the flat universe (K, = 0) and use Cartesian coordinates for simplicity. If we fix 6, 
Eqs.fl2.23j), fl2.18j), and (J2.2l|) become respectively 



1 ~ - i 4ttG , _j 8tvG . 

—DiD'a = -^-a(S\ + E) + a{a 13 a i3 + —^-pb) 



I2ttG 
— —f> b 



Id. 



2a- -0 



1 d - 
c at 



2 1 - 2V> 



1 , , AmGa 2 „, ,. 2 (1 - 2ib) 2 a 2 

lP,ki:, k = j— E(l - 2lA) 2 + - l 



-Arh 



frtij 



1 ~ 6+ 



4 

DiDjOt 



189 

c dt c a 

16ttG 



-Pb 



4 a 
1 

c a 



2 2 



-jijD k D k a 



-(q + 2)- 



+ a 



( 6 ~ 



8tvG / £ 



where 



R 



1-2-0 
1 



(V>,y +SijA f ip) + 



(1-20)' 



1 Q i 



0(hip,h 2 



20 dyK «'' 



(2.30) 
(2.31) 

(2.32) 

(2.33) 
(2.34) 



and Af is the Laplacian with respect to 5ij. 

Finally, we give the equations for matter. Since we shall restrict ourselves to the evolution of the density fluctuation after 
the decoupling of matter and radiation, we shall adopt pressure-free dust as a model of the matter, which will be a good 
approximation to the present universe. The energy momentum tensor for the dust is written as 



Tn,L> 2 U. V 

= pC U U , 

where u M and p axe the four velocity and the density, respectively. The p obeys the continuity equation 

= o. 

Also the relation between u 1 and u° becomes 



dx 1 
~dt 



dx'/dTp u l 



(2.35) 
(2.36) 

(2.37) 



cdt/drp 

where t p represents the proper time of a dust element. We should note that v l represents a peculiar velocity of a fluid element 
in comoving frame, so the physical peculiar velocity is av l . Using the peculiar velocity v 1 , the explicit form of the continuity 
equation becomes 

1 dp, , ld{p,v l ) 
c~df 



+ -- 

c 



= 0, 



(2.38) 



where p* = ay 2 a 3 pu° is the so-called conserved density ( Chandrasekharl965 ; Willl992 ) . The equations of motion are derived 
from 



T M , ;M = 0. 

Making use of Eq.j2.38j), it takes the following form which will be convenient for our present purpose: 



1 dui v 3 dui 
c dt c dx 3 



-aa,i u° + f3 3 ,i Uj 



2a 



1_-jk UjUk 
7 ..n ' 



(2.39) 



(2.40) 
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2.2 Cosmological post-Newtonian approximation 

There have been many studies of the post-Newtonian approximation in the cosmological circumstances in the Eulerian 
framework ( [Nariai fc Uenol960| [trvinel965| ; |Peeblesl98(j| ; |Tomital984 |Futamasel9"8^; |Futamasel989|; [Tomital99l[). Here w e 
shall develop the PN approximation based on the (3+1) formalism described above ( Futamasel99(j ; Shibata & Asadal99E 



We remark that the metric required to the usual Eulerian Newtonian picture in the perturbed FRW universe is known 
to take the following form: 

ds 2 = - (l + ^) c 2 dt 2 + a 2 (t) (l - ^) 8 tJ dx'dxi, (2.41) 
where (f>M is the Newtonian gravitational potential related to the matter density fluctuation field 5m via the Poisson's equation: 



A f <t> N (x,t) =4TiGa 2 (t)p b (t)S N (x,t) 
where pb = {p(x,t)) v and 



(2.42) 



5 N (x,t) 



p N (x,t) - p b {t) 
Pb(t) 



The above metric is usually used to give an accurate description of the trajectories of nonrelativistic fluid elements such as 
CDM inside scales much smaller than the Hubble distance off -1 (Newtonian limit). 

We shall expand the basic equations by using c _1 as the perturbation parameter in order to have the post-Newtonian 



approximation under the condition that the lowest order of metric takes the above Newtonian form (2.41). Thus we adopt 
the coordinate conditions where the lowest order of a and ip are of the order c~ 2 , and other metric perturbations are of the 
order c~ 3 at most. Therefore, according to Chandrasekhar's description ( [1965| ), the PN terms of all metric variables used in 
this paper may be expanded as follows: 

« = l + ^ + ^ + 0(c-% 

,/,(2) ,14) 







+ 0(c" 



h (4) 



+ 0(c~ 



f (3) 
'J 



+ o( c - 



0(3) 



0(c 



(2.43) 



where <f> is the Newtonian gravitational potential as shown below. Note that t/)*- 2 ' is not same as a priori, but actually we 



will see that they coincide. Thus the metric up to the order c~ agrees with ( [2.4l| ) . Also if we assume that the lowest order in 
P' is of the order c -3 , the lowest order in and 9 both become c -3 through Einstein equations. Using the above variables, 
the four velocity is expanded as 



u 1 I 

u — 1 + 



u 



{.i(!« v -*)} +0 < 

, + {i(i«v + *)} 
1+ {i(l«v-,)} 

«* f w* a 2 : 

h \ -r I —ol v 

c \ c 3 \2 



+ 0(c- 



0(c" 



- 2i, <2) ) + ^-p 



1 



+ 0(c" 5 ), 



(2.44) 



where terms in the bracket {} denote the first PN corrections, and v 1 is equal to that defined in Eq.(2.37) and v 



Following Shibata and Asada (1995), we consider only the PN expansion of /3 1 , not of since only f3 l appears in the above 
expressions. E,Ji,Sij and 5*'; are likewise expanded as 



E = pc 

7 2 2 

Ji — pc a 

a 2 4 

bij = pc a 



l + ±a 2 v 2 + 0(c- 4 ) 



V l+ V - ( a 2^2 _ _ 2^(2)") + ^^(3), + 0(c - 5) 
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2 2 

pc a 



1 2 . -4\ 
-zV +0(C 

c 



(2.45) 



By substituting the above expressions into Einstein equations, we can derive the relations between the metric variables 



and the matter variables in each order of c n . From Eqs.(|2.30|) and (2.45), we find 



:A f <t>+ — (A /Q (4) + 2^ 2) A f <j,-^ 2 \ k <f>, k ) 



r 47rG(p - p b ) 



1 

+ -r 



:-' , 9 n(S) ■" ■ 



4irG(p</> + 2p b( f> + 2paV) - —6 W> - 2- 

at a 



(3) 



+ 0(c- 6 ). 



The equation for <f> is derived in the c 2 order part of the above equation: 
A f 4> = 4mGa{p- p b ). 



(2.46) 



(2.47) 



Thus we may call <f> the Newtonian- like gravitational potential. Similarly from Eqs.(2.31) and (2.45), we find the equation for 
i// 2 - 1 as follows: 



A/V (2) = 4irGa 2 (p- p b ). 
So we can conclude 
i// 2) = <j>. 



The PN potentials and ip^ are obtained by taking the c 4 part of Eqs.(2.46) and (2.31), respectively, 



A /Q (4) 



</>,k +4nGa 2 (2pa 2 v 2 - p(f> + 4pb<j>) ~ °^ ~Q~l 



89^ 



2aa6 



C3) 



A/V> (4) = 4-irGa 2 [pa 2 v 2 - 4(p - pb)<t>) - ^4>,i <f>,i +c 



o(3) 



where we used the Newtonian order equations (2.47) and (2.4£). 



The relation between f3^ t and the matter variables is obtained from Eqs.(2.19) and (2.45) to be 



- Ft 1 -' 

1 ' Jl * V dt 



+ ■ 



l&nGapv 1 



(2.48) 
(2.49) 

(2.50) 
(2.51) 

(2.52) 



where we used the following relation derived from (2.20): 



(2.53) 



Finally, we give the material equations up to the first PN order. From Eqs.(2.38) and (2.40), the continuity equation and 
equations of motion become respectively 



d_ 

dt 
d_ 

dt 



^ a3 { 1 + ^(^ a2 " 2 - 3 ^)} 

'{l + i(i„V-3*)} 



+ 



_d_ 

dx i 



2 

a v 



+ v J 



pa 

d 



^{l + l(iaV-3,)} 
aV{l + i(iaV-3*)} 



+ O(c" 4 )=0, 



(4) 3 2 2, 2 a (3)j j i 9 I 2 a (3)n , i 9 . 2 a (3)i\ 

a y >,i +-a v <j>,i -a f3 y > J ti v J + — (a /3 y 1 ) + v J — - (a f3 y ' ) 



+ 0(c~ 



(2.54) 



(2.55) 



The lowest order in these two equations reduces to the Newtonian Eulerian equations of hydrodynamics on the FRW back- 
ground for a pressureless fluid ( Peebles 1980] ), and terms of the order c~ 2 provide the first PN corrections. For convenience, 
Eq.(2.55) can be rewritten as 



dv l i dv % a j 1 i 



dt 



dxi a 



&{2 aV ~ 3 V +V dxl\2 a 



1 1 



2 v 2 - 2,<f> 
d 



* (4) „ +3#,< +a 2 v 2 ^ + ^-(a 2 ^) + oVG8 (3)1 « -I3 (3)j ,i ) 
dt 



+ 0(c~ 



(2.56) 



Eqs.(2.54) and (2.56) are our basic equations for the relativistic Lagrangian approach to the trajectory field of matter fluid 
elements. 

We note that the above equations do not contain the A term explicitly, except in the equation to determine the expansion 
law, ( |2,29| ). Therefore, they may be used for the case A / as well. 

Before going into the Lagrangian perturbation theory, we consider linearized theory in the Eulerian picture for PN 
approximation on the Einstein-de Sitter background (/C = A = 0). In doing so we remark the order of linear terms of <^>/c 2 
and a' 4 '/c 4 , respectively: 



<f> Gp b (al) 2 
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*<*> G Pb (alf 4> 



(2.57) 



where I denotes the characteristic length of a density fluctuation in the comoving frame and we have used Eqs.(2.47) and 
(2.50). This equation (2.57) means that c/ 4 ' may become larger than the quadratic terms of perturbation quantities such 
as 5 2 for fluctuations with scales equivalent to or beyond the horizon scale ct. Therefore, we will leave all linear terms of 
perturbation quantities in th e abov e equ ations regardless of the power of c. 
By linearizing Eqs.( 2.54 ) and ( [2.56| ), we get the following equations: 



« («-f H-o. 

— v + 2-v = --zP.i ~ 
Ot a a 2 



1 1 



J 4 ) 



+ J(a 2 /3 (3)I ) 



where we used the fact that pb(t)a(t) is constant. Using Eqs.(2.47), (2.5C) and (2.53), Eq.(2.59) becomes 



— v ,i + 2-v , 
Ot a 



-AnGp b S - \4nGp b 2(/> - 



dt 2 a dt 



(2.58) 
(2.59) 

(2.60) 



By substituting Eq.(2.58) into Eq.(2.6C), we can get the evolution equation for the density contrast in our coordinate system: 

d 2 s 



— + 2— 

dt 2 adt 



AirGpbS ■ 



\4nGpt ■ 2(f> + 3 " 



a dt 



(2.61) 



As we will show in the next section, (j> remains constant during a period when the linear approximation is valid. Thus we can 
integrate Eq.(2.61) to get 

2/3 



/ I \ 2/3 / 2 \ 2 



(2.62) 



where ti denotes the initial time and we have considered only the growing mode because the decaying mode (oc t~ ) plays 
physically no important role. It is remarked that (S + 2<f>/c?) is a gauge invariant quantity in the PN approximation and 
is guaranteed to be a physical density fluctuation (see Appendix A). In next section we will use Eq.(2.62) as the bridge 
connecting between the Eulerian perturbation picture and the Lagrangian perturbation picture in order to determine the 
initial conditions of the first PN displacement vector. 



3 LAGRANGIAN PERTURBATION APPROACH TO THE PERTURBED FRW UNIVERSE IN THE 
FIRST POST-NEWTONIAN APPROXIMATION 

3.1 Post-Newtonian Lagrangian evolution equations for the trajectory field 

In this section our purpose is to construct the complete set of evolution equations for the trajectory field of the the fluid 
elements in the post-Newtonian approximation by extending the formalism in Newtonian theory developed by Buchert. 

In the following discussion, we only consider gravitational instability on the Einstein-de Sitter background (K, = A = 0) 
for simplicity: 

H 2 (t) = ^nG Pb (t). (3.1) 

The generalization to more general FRW background is straightforward and will be presented elsewhere. 

In the Lagrangian description, we concentrate on the integral curves x = f(X,t) of the velocity field v(x,t): 



1L = 2L 

dt dt 



= v(f,t), f(X,t T )=X, (3.2) 



where X denote the Lagrangian coordinates which label fluid elements, x are the positions of these elements in Eulerian 
space at time t, and ti is the initial time when Lagrangian coordinates are defined. It should be noted that we introduce the 
Lagrangian coordinates on the comoving coordinates because we have already derived basic equations by using FRW metric 
defined with the comoving coordinates. As long as the mapping / is invertible, we can give the inverse of the deformation 
tensor fuj which is written in terms of variables (X,t): 



OX 1 



= hi,j(X,t) — —(.iab£jcdfc\afd\b, (3-3) 



where J is the determinant of the deformation tensor f t \j, and tijk is the totally antisymmetric quantity with 6123 = +1. The 
comma and the vertical slash in the subscript denote partial differentiation with respect to the Eulerian coordinates and the 
Lagrangian coordinates, respectively. The following equations are naturally satisfied: 
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hi, k fi 



k J 



Si j, 



fi\khk,. 



(3.4) 



where Sij denotes Kronecker's delta function. Using the above equations and the following identities 



d_ = d_ 

dt ~ dt 
1 



J{X,t) dt 



d_ 

- ~ dt 
J{X,t) 



+ v 



_d_ 

dx 

dX k 



dx 1 



dv\X,t) 



dX k 



the continuity equation ( 2.54 ) in the PN approximation may be written as 



pa 



l+-A(X,t))J(X,t) 



d_ 

dt 
where 



A(X,t):=^( t) (%) 



0(c 



= 0, 



(3.5) 



(3.6) 



(3.7) 



3<j>(X,t). 

Thus the density field is integrated exactly up to the first PN order in the Lagrangian picture as in the Newtonian case: 
p(X, t)a z J{X,t) = (l + ^A(X,t)) ~* P (X) a 3 (l + 1 A (X)) + 0(c- 4 ), (3.8) 
where the quantities with ° such as a are the quantities at the initial time ij, and 
A (X) = ±°atf (X)-3°4>(X), 

A x I (X) = 4ttG a (p (X)-p b ), 

with Ax being the Laplacian with respect to the Lagrangian coordinates. We should remark that the p in all equations below 
will appear in terms of J p. 

Next we consider the trajectory field / of the pressureless fluid elements in PN approximation. The field / is given by 
solving Eq.(2.56). The equations of motion can be rewritten formally as 



d 2 P n adf 
— - h 2- — 

dt 2 a dt 



< N V(X,t) + ^ PN) g{X,t) + ^ PN V(X,t) + • 



(3.9) 



where ^ 'g % , ( PN ^g l and ( 2PAr 'g I denote the acceleration fields in the order c° ,c -2 , and the order c -4 , respectively. By means 
of the perturbation theory, the trajectory field given by integrating the above equation may be written formally as 



f(X,t) = X 



^p(X,t) + ±^ P (X,t)- 



(2PJV), 



p(x, t) ■ 



where p denotes the displacement vector and the following initial conditions must be imposed by Eq.(3.2): 

<%(x,t,) = (p %(JM/) = -.. = o. 



(3.10) 
(3.11) 



Now we construct the set of equations to solve the trajectory field /. Likewise as a Newtonian case ( Buchertl995 ), we 
consider the divergence and the rotation of equations of motion in PN approximation with respect to the Eulerian coordinates, 
respectively. In the Lagrangian picture, the former set becomes basically the evolution equation of the lo ngitu dinal part of 
the trajectory field and the latter becomes one of its transverse part. First by operating divergence on Eq.( 2.56 ) with respect 
to the Eulerian coordinates and using Einstein equations in each order of c, we obtain 

\d 2 p 



d 

dx l 



dt 2 



a dx 1 



dt 



« P (i + ii)f«»-Ii- 

!♦)}■ 



' df l dA~ 


1 ■ 


~dt~dt 





+[v 



H f 2 (d±, 

a 2 dt I \dt 



( d<t>,. 



+ v 



(3)i 



(3)j 



+0(c~ 



(3.12) 



According to the usual procedure in the Lagrangian formalism, we modify the above equation by using Eqs.(3.2), (3.S) and 
(1.4) in the following form: 



2 €jab£icdfc\afd\b 



a dp 
-)- 2 — 

dt 2 a dt 



„ 1 o o 3 

-4ttG^t Pa 
a s 



A) +4nGp b J- 



d 



2c 



2 £jab£icdfc\afd\b q 



af_dA 
~dt~dt 



c 2 J 



€jab£kpqfc\a fd \bfc\p fd | < 



dfi\j dfi 



2^ + ^f^k 



4ttG 



'dfi 
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3J d 
a 2 c 2 dt 
J dfi ( d 



/ 2 fd4> 1 , . , dfi a \\ 



Jdfifd 1 d/ fc 8 d\ 1 ( . 

+ C 2 dt \dt 2J £labekcdMaMb dt dX' + a) 2 J €mpqeirsIr}pMq(plr ' 



(/? (3)i ,i - /3 (3)J , 1 ) +0( C " 4 ), (3.13) 



where 



r ] , / v- . (JV) . 1 (PN) 

J = det I a + 'p + — K 'p + 



and we have used Eq.(3.8), which was obtained by solving the continuity equation, and the Eulerian partial time derivative 



of 4> was rewritten as, for example, 



dt dt dx 

d<j> dfi 1 



dt dt 2 J^ ab€ ' cdf ^ afd ^- 



We did not transform the terms including the shift vector f3^ 1 in Eq.(3.13) into the form in Lagrangian picture on purpose, 
because its expression in Eulerian picture will make us clearer to understand the expansion of the above equation introduced 
in later discussions. 



Next by considering the rotation of Eq.(2.56) with respect to the Eulerian coordinates, we get the following equation: 

d 2 f h 



d 



dt 2 
1 1 



2 dd£_ + 1 df dA 

a dt c 2 dt dt 







d 



2r,(3)k\ . 2 i 



0{c 



(3.14) 



Similarly we again modify the above equation by using Eqs.(|J), Q and Q in the following form: 

€abcfj\afi\b 



d fj\c r,d d fj\c 

dt 2 a dt 

1^/2 a (3)k 



.2 dfi dfi\ c 



(3)k 



-^£abcfj\afi\b 



,i — 
dA 



P (3)l ,k)} 
dfj dA\ c 



+ 



+ 0(c" 4 ), 



dt dt dt dt 

where the term containing the Eulerian partial spatial derivative of (j> was rewritten as 
9 i 2 , \ 2 df m df m \ c 



(3.15) 



dt dt 



We note that the right-hand side of Eq.(3.15) is of the order c -2 , namely, has only the first PN order quantities since we 
consider the gravitational instability of the self-gravitating dust model. Eqs.(3.13) and ( 3.15 ) are our basic equations for the 
trajectory field /, in the Lagrangian picture in the PN approximation. These are evolution equations for general trajectory 



fields which can be applied even in the nonlinear situation (S ^> 1) as long as the PN approximation introduced in §2.2 is 
valid. The structure of these equations allows us to solve the displacement vector p iteratively up to any desired order in c~" 



in terms of the lower order displacement vectors and metric variables. (The shift vector is expressed by solving Eqs.( 2.52) and 
(2.53) in terms of Newtonian order quantities such as <j>.) Thus we need to know only the solution at the Newtonian order to 
have the solution in the first PN approximation. 



3.2 The irrotational flow for the trajectory field in PN approximation 

In this section, we investigate the irrotational flow for the trajectory field. 



3.2.1 The basic equation 



When we consider irrotational flow up to PN order, we need only Eq.(3.13) among the above basic equations. Further 
simplification is available in this situation. In Newtonian case, only an irrotational solution survives among general solutions 
by employing the appropriate initial conditions that the initial peculiar velocity field is proportional to the initial peculiar 
gravitational field on basis of the results derived by Newtonian linearized theory (Buchertl993; Buchertl994). Thus we consider 
the following trajectory field with the longitudinal perturbations: 

1 



f(X,t) = X 1 + (N) S lt (X,t) + 



(PN) 



'S ll (X,t) + ^ 



(3.16) 



Note that the following initial conditions must be imposed: 
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(N) 



, 5| 1 (X,t I ) = ( pjv ^| l (X,t / ) 



0. 



It should be noted that the above expansion of the trajectory field here is a post-Newtonian expansion introduced through 
the equations of motion in PN approximation as discussed in the previous section, and we impose no restriction on the 
amplitude of density contrast field. Neglecting the vorticity is known to be a reasonable assumption for the pressureless fluid 
in Newtonian theory, but this is not the case for the PN displacement vector. We shall discuss the PN transverse flow in our 
coordinates later and here concentrate only on the irrotational flow. Inserting Eq.(3.16) into Eq.(3.13), we obtain the following 
equation from the c° order: 



"^tjab^icd Jc\a Jd\ 



2 

where 



47TG o3 /o 

— a V 



Pb Jn 



(3.17) 



(JV) 



%(X,t)=X^ + ^%(X,t), 



Jn = det 



det (Sa + (N) S 



and the dot denote the Lagrangian time-derivative. In this case, Eq.(3.1'i) entirely agrees with the result derived from the 
following familiar set of equations in Newtonian theory in the Lagrangian description, 



Pnci 

£ 

dt 

1 



a 



o o 3 

- ?!L 

Jn 



(It 

(N). 



2J, 



-e ja be,J N \f c J N \f, 



d\bt>N\jy 



i \ i ,.. (JV) r 9 
^Cjabtkpq Jc\a Jd\b Qj£j 



( 1 Wf , (JV) 

\Jn Mp ' 



d\q<PN\k 



\ 4ttG /o 3 3 

J = — — I Pa -J N p b a 



that is, these equations can be changed to the following system of equations in the Eulerian description: 
d{pNa A ) + d(p N a 3 v N ) = Q 



dt 



dv N 



~KT + v N-^-T + 2 ~ V N = 5^ 



dxi 



A x (j} N = AnGa (p N - p b ) 



(3.18) 
(3.19) 
(3.20) 

(3.21) 

(3.22) 
(3.23) 



where it is noted that a physical peculiar velocity of a fluid element is a{i)v l N . 

Now, we derive the equation which governs the first PN irrotational velocity potential ^ PN ^S. For this purpose we must 
express terms in the right-hand side of Eq.( [3.13[ ) explicitly by Newtonian quantities. This will be easily done by using Eq.(3.£ ) 
and assuming 



Jn > -ttJpn, 

where Jpjv is a term of the order c -2 in expanding J: 

? , _ (PN) q , ((PN) q (JV)o (PN) q (JV)e \ , 1 (PN) q (N) q (N) q 

JPN — ~T y >J\ii ~ J I T^^iab^jcd d\ca *->\db »-> 



(3.24) 



(3.25) 



(The assumption ( 3.24 ) seems to be satisfied almost everywhere in the situation we are interested in here, but it may break 
down at a shell-crossing point, where the Lagrangian perturbation theory itself also do. We leave this problem for the pr esent .) 
According to the above discussion, for example, we shall make the following replacement in the right hand side of Eq.( |3.13j ): 

1 , , f f f , dfm dfi 

2 7 3 a b^kpqjc\ajd\bjc\pjd\q *^ dt 



c 2 Jn 
1 

2c 2 J N 
1 



tjab^kpq 



J c\a J d 



(JV) f (*), f^i rf(JV) /' 
b Jc | p Jd | q 



dt 



dt 



\ k +0(c- 4 ) 



(N) p (N) r {1\ ) r {1\ ) r 
tjabticdtirstkpq Jc\a Jd\b Jr\p Js\q 



(N), 



w f (N )f f%i d(N) f' d 

c2 tjabe lcd J c]a f d \b dt dt dt 



dt 

2 d (jv) /r 



dt 



dt 



0{c 



!>N\k 



0(c 



(3.26) 



where we have used Eq.(3.19). Similar calculations lead us to the equation for the first PN displacement vector with the 
longitudinal perturbation we are looking for: 



1 (JV), (JV)r 

^tjabticd Jc\a Jd\b 



(PN) q ,r, a (PN) q 



AnG o o 3 o 

5- Pa A +ATvGp b 



Jp N + a 2 ^ N %yj N 



+ II( iN) S, (piv) S) + M(^S, W S, (PN) S) 
3 



47rG o 3 

— s- Pa 



5<f>N 
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+3Jjv 



1 



i <j>N a #jy 
dt 2 a dt 



+Jn 
' 9_ 



d(a 2< - N % 

dt 



' (N) . dVW%) d (iv) . rf(a 2(Ar) 5|,) 



( W S| 1 i,)+( (, "S |fc 



2 d^^S,, 

) ir 



(JV)A (JV)c 



eft 



dt 



(3.27) 



where 

{PN) S) = (PAr) 5| 



2- (iv) S, 
a 



\ _ (PiV) c f (N) 



(N) 



'S Hi + 2 



lu 
a 



"t~ 2— ♦S'jij 



'9, 



(3.28) 



Here, as expected, all quadratic terms on the right-hand side are expressed only by Newtonian quantities. Since we have so 
far assumed that Newtonian flow is irrotational, the f3 l N terms in Eq.(J.13) are canceled out in Eq.(3.27), that is, we can 
conclude that the gauge conditions do not have an influence on the trajectory field in the first PN approximation. On the 
contrary, the assumption (3.16) may mean that we tacitly impose gauge conditions that there exist coordinates which have 
the hypersurface normal to a congruence of geodesies of the pressureless fluid elements. Eqs.(3.17) and (3.27) are the master 
equations for solving the longitudinal perturbations ' 'S and (PN ' I S, respectively. It should be again noted that since the 
continuity equation is exactly solved, these equations can be integrated for some initial conditions only if we give the initial 
density field P as the source function. 



3.2.2 Longitudinal perturbation solutions up to 1st PN order 



Let us now solve the equations Eqs.(3.17) and (3.27). We first consider Eq.(3.17) which is a master equation for solving ^S: 



AN) 4 



— tjabticd J c \a Jd\b 



4ttG o3 / 

— a I 



P - Pb Jn 



The Einstein-de Sitter background solution of Eq.(3.1) can be chosen as 

2/3 

a(t) 



( L\ 



(3.29) 



(3.30) 



where to denotes the present time. We use this normalization for the scale factor since it makes the physical interpretation of 
the solution more transparent partly because the length in the comoving frame is equal to the physical length at the present 
universe with this choice. We shall assume that the initial hypersurface is sufficiently early so that the initial density contrast 
field is much smaller than unity. This allows us to adopt the order of the density contrast, say, A as a new perturbation 
parameter. Thus we calculate the following ansatz for the longitudinal perturbations superposed on an isotropic homogeneous 
deformation: 



'fi =X % + C'Sij = X' + Xq z (tM>(X) + X'q zz (t)^>(X) + 



„(2), 



(3.31) 



where the fist term on the right-hand side represents an isotropic homogeneous flow for the FRW background. The parameter 
A is assumed to be small and dimensionless. As described above, it can be considered as the amplitude of the initial density 
fluctuation field. We formally split the initial density field accordingly as 



P(X) =p b +X P b S(X) 



P(X] 



(3.32) 



where $ denotes the initial density contrast field. This choice of the initial data is adequate, since the density need not 
be perturbed in the Lagrangian framework. For simplicity we make use of the usual assumption that the initial peculiar 
velocity is proportional to the initial peculiar acceleration (Buchertl992; Buchertl994). This assumption has been proven to 
be appropriate in Newtonian linearized theory (see e.g., ( Peeblesl980| )). The Newtonian peculiar velocity field' iv 'u(X, t) and 
acceleration field ( - Jv ^(X, i) is defined as the following forms in terms of the trajectory field f(X,t) on the comoving frame: 



iN) u(X,t) = a(t)^f(X,t), 

w w(X,t) = 2a(t) {Ny f{X,t)+a{t) {N) f(X,t) {= —^n) 
Therefore we employ the following assumption for the initial data 
u (X) = tiw (X), 

so the initial conditions for q z <p and q^zf become 
q^t^iX) = 0, q„Qti) V $\X) = 0, 



(3.33) 
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Xa(t I )q z (ti)^ ) (X) = (N) w(X)t I , q zz (t I )<p\i\X) = 0, 

X(2a(t 1 )q z (t I ) + a(t I )q z (t I ))^(X) = W w(X), (3.34) 
where ^ w denote the initial peculiar acceleration field constrained by 
Vx ■ (N) w= -(1 + ZI )A X °<t>= -2/(3t?(l + zi))X 8 

because of a= 1/(1 + zi) where zi denotes the redshift of the initial time. With a superposition ansatz for the Lagrangian 
perturbations of an Einstein-de Sitter background we have obtained the following family of trajectories ^ N 'f as an irrotational 
solution of Eq.(p\29|) up to the second order: 



t \ 2 / 3 



/3\ 2 r 3 / t \ 4 / 3 3 / t \ 
*" = UJ -14 U + -AT i ) 



t\ 2 /'i i 
2 



35 \ti 



and 



Ax<p 



(2) 



2 A 



Note that cj> is the following order: 



(X) = —A 



(i) 



t?(l + *i) 



(3.35) 

(3.36) 
(3.37) 

(3.38) 



Let us now turn to the longitudinal solution for the PN displacement vector. First of all, we need to express the Newtonian 
potential 4>n in terms of the known trajectory field ^ N 'f because Eq.(3.27) for solving the PN displacement vector includes 
(j>N in the source terms. Using Eq.(p.3|), Eq.(3.19) can be rewritten as 



(N) 



'hj,i(f>N\j = 



}d (N) f i 



(It 



and multiplying ^ N 'fi\k which is inverse of Whi i , we get 

2d (N) f i 



, _ (]V) , - 
N\k — Ji\k — 



d_ 

dt 



dt 



(3.39) 



Moreover, by differentiating the both sides of Eq.( p.39| ), we can get the Poisson's equation with respect to the Lagrangian 
coordinates: 



Ax4>n = —A (2aaq z + a 2 q z ^j Axtp W — A 2 



(2aaq zz + a 2 q zz ) A x ^ 2) + q* {2aaq z + a 2 q z ) (v^jf^J 



+ 0(A 3 ) (3.40) 



Therefore we can arrive at the expression of 0jv with respect to the Lagrangian coordinates: 



(X,t) 



-A (2aaq z + a q z ) ip 
-X 1 



' (X) + q z (2aaq z + a q z ) 4>(X)] + 0(A 3 



where 



(X) 



dY 



*?(l + z/) 2 
4tt|X -Y\ 



m (X)-X 2 [(2adq zz + a 2 q : ^ (2) 
^(i) _ ^2 ^2aaq Z z + a 2 q zz ^j i/? (2) + q z (2aaq z + a 2 q z ^ (j>{X)\ + 0(A 3 ), 



(3.41) 



(3.42) 



Here it is easy to check that the Newtonian potential 4>N satisfies the initial condition Eq.(3.38). It should be noted that the 
order A part of (j>N becomes constant because of considering only the growing mode of q z , which is consistent with the result 
of linear theory on an Einstein-de Sitter background as discussed in §2.2. 

Now to get the evolution equation for ' p 'Su, we consider Eq.(3.27). Similarly we make the following longitudinal 
perturbation ansatz: 

(PN) S {1 = XQ z (t)^ ) (X) + X 2 Q zz (t)$\ 2) (X)+0(X 3 ). (3.43) 

We m ust consider the initial conditions of Q z & x \ Q zz & 2 " 1 and so on in order to solve the second rank differential equation 
(3.27) according to the above ansatz. Part of the initial condition must be imposed by the definition of the displacement 
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vectors: 



Q z (t I )$W(X) = Q, Q zz (t I )$\ 2) (X)=Q,--- (3.44) 

Although we need another initial condition for the peculiar velocity field, we leave Qz{ti)& and Qzz(ti)& and so on 
arbitrary for the present and think about them later. 



The order A 1 part of Eq.(3.27) becomes 

2 2ip w (X) 



a 3t 2 



A X $ (1, (X) = A 



3t 2 (1 + Zl )H 2 ' 



Hence we have 



Qz+2-Qz 
a 

with 

Ax* (1) (X) 



— Q =- 

3t 2 ^ z ~ t 2 



(3.45) 
(3.46) 
(3.47) 



It may b e noted here that defined by Eq.( 3.47 ) is the same as so-c alled "superpotential" on the hypersurface at the initial 
time ti ( ChandrasekharI965 ) . We can find the general solution of Eq.( 3.4(: ): 



,2/3 3 i 

at ' - - + c 2 t , 



(3.48) 



where ci and C2 are constants of integration. To determine them, we make use of the density filed (3.8) solved up to the first 
PN order, and consider the linear stage where qip^ and Q&^/c 2 are much smaller than unity. Then the linear order in A of 
the density field (3.8) becomes 

5(X,t)=°5(X)-qz(t)A xlf w (X)~ ±-Qz(t)A x <S> (1) (X) + 0(\ 2 ), for q z (t)A xV [1 \ ^Qz (t)A x $ w < I, (3.49) 



whe re we used the fact that <j> =<f> in the order A 1 through Eq.( 3,41 ). By substituting Eqs.( 3.35 ), ( 3.36 ), ( 3.47 ) and ( 3.4£ ) into 
Eq.( 3.49 ), we obtain 



S(X,t) = (±) 5 (X) + A % (X) ( Cl t 2 / 3 - | + c 2 t 



(3.50) 



By comparing this with the density contrast field ( 2.62 ) derived in the linearized theory within the framework of our PN 
formalism, we can arrive at the solution of Q z - 

t \ 2 / 3 



«.(*)= *.(*) = (§) (£) ; -i 



(3.51) 



It should be again remarked that although the density contrast is supposed to be much smaller than unity in linearized theory, 
we do not have to assume such smallness for the density contrast in Lagrangian perturbation theory. 

The initial conditions for the quantities of higher order in A has no corresponding counterparts in linearized theory and 
thus we are not able to know the initial conditions for such quantities. We may, however, expect them to be much smaller 

( 2) 

than the order A quantities at initial time so that we will neglect them and we adopt the following initial condition for 
and so on: 

Q„(ti)$\?(X) = 0,---. (3.52) 

We are now in the position to be able to solve Q Z z& 2 \X). The differential equation for Q Z z^ 2 ^ becomes the following form 
by using the fact that 2aaq z + a 2 q z = + Zi) 2 tf) 

2 , 



+ 2-Q z 
a 



+ 



3t 2 



A X & 2) (X) 



(Qz + 2^Q Z ) - Qz (q\ + 2-^ + ^i-QzQz 



.a „„... „ da 

> h 12aa + 2a — - 

a at 6 



(1) i 1 1) 

IP, . . — 9). .CO,. 



■2 \ ■ • d 3 q Z z 2 



qzz + 6 (3d + ad') q zz + 9ad- 



3 \ q z + 3^q z J (2adg 2 + a 2 g' z ) + ^<2z (2adg z + a 2 <j' z ) 
V W (X 



I- + a 2 q z y ''— ^g 2 {2aaq z + a 2 q z ) + 



3t 2 



dt 3 
X) 
1 



+ a' 



(X) 

d 4 q Z z 
dt* 



[2aaq Z z +a qzz) 



v {2 \x) 



2{1 + ZI )H 2 2 



5 2-2 



+ w 2s{x '(l+ ZI ) 2 t 2 



(3.53) 
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By using Eqs.(3.35) and (3.51), this equation is rewritten in the following simple form: 

S(2), 



+ 2-Q 

a 



3 / M~ 2/3 
5 \ti) 

2 



Ti) ~ 5 



+ 



3tf (*J 



(1 + Zl) 2 * 2 



2*(X)^(X)-I(^(X)) 



2 -^(X) + |^(X) 



(3.54) 



To obtain solutions of Eq.(3.54), we make use of the linearity of Poisson's equation and split the potential $^ 2 'into two parts. 

(3.55) 



Then, we have to solve separately the linear ordinary differential equation: 



)(<•) , o£n<» _ A r>(") _ _L 

tzz +^ a ^zz M2 ^zz 



3 / t \ -' 2/3 ( t \ ~* 2 / * \ 



2 / t \- 7 ^ 



with 

A x 4> (2a) = 2 
and 

q£ + 2^?2 

a 

with 



Ax$ (2b) 



)(<0 



3i 2 



-f-r 2 

3* 2 UJ ' 



(1 + z,) 2 * 2 



2 <S 



15, 9 



(2) 



We find the following solution of Q zz for the initial conditions Eqs.(3.44) and (3.52) 



Q 



(b) 



(2) 5(4/) 



2/i\ 2 / 3 2 £ 
3 + 15 



3 / t \ 4 / 3 3 / t \ 2 / 3 1 4 / t \ - 1 
li\Ti) + 5V*7/ ^ 2 + 35 Vt7/ 

(5)"' 



(3.56) 



(3.57) 



(3.58) 



(3.59) 



(3.60) 



We have obtained the following family of trajectories x = f(X,t) as an irrotational solution up to the second order of A in 
the perturbed Einstein-de Sitter background: 

f(X,t) = X + q z (t)V x (<p W (X) + ^<& (1) (X)) +g zz (t)Vx(V 2) (X) + ^<& (2a) (X)) + ^Q z b) z {t)V x ^ 2b) {X), (3.61a) 



with 



(!) (£) 



£ \ 2 / 3 



.3,2 

(!) 



/3V 3 / i \ 4 / 3 3 / * \ 2 / 3 1 4 / * \ 
^ = Uj "14 1*7 J + 5UJ - 2 + 35 1*7 J 



5\ti) ~ 3 + 15 



and 



Ax<p 



(i) 



2 o 
"3* 



3(1 + 2,) 2 * 2 



Ax$ (2a) =2f$( 1 V[ 1) -$« 1) ^ 1) 
l 



(2b) 



(1 + z,) 2 * 2 



2,5^ 



(i) 



15 7 9 



(2) 



(3.61b) 
(3.61c) 
(3.61d) 

(3.61e) 
(3.61f) 
(3.61g) 

(3.61h) 
(3.61i) 



3.2.3 Physical Interpretation of the longitudinal perturbation solution 

In this subsection, we will express the trajectory field ( fc,61a ) in terms of physical quantities. Our coordi nate co nditions ( 2.2( ) 
show that the density contrast field 8 has gauge ambiguities. In order to interpret the trajectory field ( 3.61a ) physically, we 
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should express f l (X,t) in terms of the gauge invariant quantities at the initial time, rather than the coordinate dependent 

o 

initial density field 8, because only gauge invariant quantities are physically relevant. This will turn out to be essential for the 
physical interpretation. The physical density contrast becomes at the initial time ti 



2 o (2) o 



8 g (X) =8+- 
where 



A x a (2 \x) = A7vG apJ g (X) 



5 {X) + -4>{X) 



(3.62) 



and S g denotes the gauge invariant quantity at the initial time defined by linearized theory, and a 1 - 2 ^ /c 2 = a — f. ft is supposed 
that a' 2 ' remains constant with respect to time in the linear regime because of neglecting the decaying mode. (The reason for 
this choice is described in Appendix A.) Thus, the quantity 8 has the spurious mode 2 (j> jc caused by our gauge conditions 
(2.26), so we had better redefine the displacement vector of ( 3.61a| ) in terms of S g , that is, since we have 



A x (^ + ^\X))=- 2 -(8(X) + lhx. 



A x ^ = ~ 8 g (X). 



(3.63) 



we may redefine the first order velocity potential in Newtonian theory as 

2 
3 

Naturally, the quantity * at is a gauge invariant variable. Thus the lowest PN irrotational velocity potential in our 

o 

coordinates is a spurious mode generated by the gauge ambiguities of the density contrast 8 in our coordinates. Therefore, it 
has no physical meaning and can be understood as a part of the gauge invariant Newtonian first order quantities. Likewise 
since we have 



Ax (V 2) (X) + 1<I> (2[1) ) 



(i) (i) (i) (i) 



\ *N\ii*N\jj 



: N\ij 1 N\ij I + ^A C 

we may redefine the second order potential in Newtonian theory as 



^■X* N — w N\jj 



N\ij N\ ij ' 

Furthermore for the first PN displacement vector <&( 2b ) we have 



(3.64) 



(l + 2,) 2 i? 



(l + zj^t 2 



2 15 
+ — 



dY 



2 8 9 ^-\{^)\^J d Y 



An\X 



2^ 



(i) 



4tt|X -Y\ 



+0(c 



so we may redefine the first order potential in the fist PN approximation as 

(*W j( Y)*W(Y) 



Ay*, 



(l + zi) 2 t 2 



2 5 9 *^-2 



(i) 



4tv\X 



(3.65) 



In the end, we have obtained the following family of trajectories x — f(X,t) as a physical 1st PN irrotational solution 
up to the second order of A in the perturbed Einstein-de Sitter background: 

1 



f(X,t) =X + q z (t)V x ^'(X) + ? «(<)Vx^(X) + ^Qf> (t)Vx*Piv(X) 



with 



Qz = 



t .2/3 



(!) IG 

(IV LiL ( lY /3 lflY /3 ~- JLflY 1 

: \2) 14 UJ 5\ti) 2 + 35 U J 



Q 



(*>) 



2 ( t \ 2 / 3 2 4 / t 



(2) 5 {h) 



3 + 15 



(3.66a) 

(3.66b) 
(3.66c) 
(3.66d) 
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and 



Ax*« = -1 S B , 



(i) 

N\ij ~* N\ij' 



(i) 



2 5 s *iv 



(i) 



+ - 



15 



dY 



(^%.(Y)^%{Y) 



4tt|X 



(3.66e) 
(3.66f) 

(3.66g) 



It should be again noted that Eq.( 3.66a) is expressed only in terms of the gauge invariant quantities at initial time when 
the linearized approximation is valid. In other words, if using trajectory field (3.66a), we can accurately describe evolutions 
of all scale fluctuations in the early stage according to the gauge invariant linearized theory, and also take into account 
the relativistic PN corrections which are generated by the non-linearities due to evolutions of fluctuations in the subsequent 
structure formation (see below in detail). 



3.3 The transverse perturbation solution at 1st PN order 



We now consider the transverse part of the trajectory field at P N order. We assumed that the vorticity for the Newtonian 
displacement vector is negligible, which seems to be reasonable ( |Buchertl99^ ; |Buchertl99j |Peeblesl98(i|) . In this section, it 
will be found that the PN displacement vector has a transverse flow with the growing mode in the Lagrangian coordinates 

even in the case. 

According to the result ( 3.66a ) in the previous section, we consider the following ansatz for the transverse perturbations 
superposed on the solved irrotational trajectory field up to the first PN order: 



f(X,t) = ( N Y( X ,t) + X 2 ± [ Qi b }(t)y PNli (X) + ± 



(PN) a {l)T 



's^ 1 (x,t) + x 2<PN) sl 2)T (x y t) 



(3.67) 



where 



{N) f(X,t) =X l + Aq z (t)*« (X) + \ 2 q zz (t)V%, 

and the PN transverse displacement vector at each order is constrained by, 



CFJV) C (2)T 



= 0. 



(3.68) 



We need only Eq.(3.15) derived in §3.1 to solve the transverse part ( PAr )5" Tl Q f the PN displacement vector: 



£abcfj\afi\i 



dt 2 a dt 

' 9 (a%^, 3 -) +e ijk £-{v> 



2 dfi dfi\c 

+ jeabcji\ b dt df 



■^ £ ^M\" \-dT dF+ dF— ) +0{ - c )> 



(3.1! 



Naturally, the order c° in this equation produces Newtonian counterparts in the case that the trajectory field is introduced 
on the comoving coordinates (Buchertl992; Buchertl993; Buchertl994). Then taking into account order of A in each term of 
the right-hand side and using the assumption ( |3.24| ), the above equation becomes 



eabcfj\afi\b 



fj\c + 2 fj\c 



= -^?|(° a ^ W ^) " ^iirf 1 ^" -P {3) '*)]+0(c-\X*) , (3.70) 



where we neglected the terms like 4>nv in our approximation because they are higher than third order A as derived in the 
previous section. This equation shows that the shift vector generates the first PN transverse part. Inserting the ansatz (3.67), 
we can rewrite the left-hand side of the above equation (3.70) as 



eabcfj\afi\b 



fj\c + 2 ~/j|c 



***** { {pN) si7 + 2 ? p xr) + ( (p; xr + 2 ? p xr) + a2 



+ 2-q z e, 
a 



i kV N\il J l\k 



+ 0(c- 4 ,A 3 l 



Accordingly, as the irrotational case, we can obtain the following differentiation equation for solving the transverse part of 
the first PN displacement vector at the order c~ 2 of the both side in the above equation: 



y&U) (PJV)o(l)T 
N\jl l\k 
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+A 



= -$|(« a ^ w * J )-Ji, 6M ^[^*u-/^ ), *)]+0(A-) 



(3.71) 



The form of this equation shows that it could be solved if we expressed the shift vector /3' 3 ' 1 in terms of already-known 
Newtonian trajectory field ( N ^f. 

Therefore, we first have to s olve the shift vector wi th respect to the Lagrangian coordinates. Following the relativistic 
linearized theory ( Bardeenl882 ; Kodama fc Sasakil984 ), it will be convenient to decompose the shift vector (3^* into its 
longitudinal part and its transverse part on the Einstein-de Sitter background: 

0(8)* = 0(3) + /3 (3) lj ^(3)^ = Q (3 72) 



Eq.(3.71) shows that we need only the explicit form of its transverse part (3^ in the situation we are here interested in. 
Since the gauge condition (2.26) implies that the transverse part of shift vector is a gauge invariant quantity, any ambiguities 
caused by the gauge freedom do not remain. The perturbation quantity f3^ 1 is constrained by Einstein equation ( [2,5' 



A/A 



(3)« 



\—Q t 1" ~ < P N ' i ) + 167r ^ a P V N- 



(3.73) 



Since the right-hand side of this equation is shown to be divergenceless exactly through the continuity equation (2.54) at 
Newtonian order, the condition (|j.72j) for 0^ is al ways satisfied. We have already obtained the expression of the peculiar 



velocity field and the peculiar gravitational potential (3.2) and ( 3.41 ), respectively. Therefore, it will be adeq uate to make use 
of the equation ( 3.73 ) for our purpose. By using Eqs.(3.3) and ( 3.24 ), we rewrite the left-hand side of Eq.( 3.73 ) in terms of 
the independent variables X and t: 



dX k d 
dxi dX k 
d 



= {N) h k 



' 3 dX k 



n.i,jP T \i 



1 (JV)r (jV)r & 

n T tkabtjcd Jc\a J d\b q-c^i. 
ZJN OA 

1 



— e, e- (iv) f , <N) f , B {3) \ 

2J N Pq JrS * P J s\qHT I 



2J% 



1 (N) r (N)r r (N)r (N)r a(3)i , j (N)r (AO r a(3)i . T (N)r (N)r 

7r-^-e kab ei pq K '} c \ a y >f d \ b -J N \ k 'f c \ p y 'fd\qP T \i+2Jn 'fc\ P k fd\ q P T \l + -JN y 'jo\p Jd\ q P T \i 



J 3 



Ax/?<? )1 - + 3Ag,*«. .4 3)i | M - 2X qz ^\ jk pP\ jk + 0(A 3 ) 



(3.74) 



where we used the fact that the perturbation quantity p\, is of order A 1 at most and 
J i v = l + A 9z *W i +0(A 2 ). 



Using Eqs.(3.8), ( p.l9[ ) and (3.66a), the right-hand side of Eq.(3.73) can also be rewritten as 



(r.h.s) 



" N dxi 



0N,i H 4>N,i 

a 



16nG o Q 3 



pa J N v l N 



d_ _ (jvj ; _d_ 
dt Jj dxi 



167rG o 3 (n] ■ 

Pa 1 >fi 



dt { 



+ 5odW/i + 2a 2 ( | + | j W/ 4 - -Le kab e 3 J N %J N % b (a^ N % k + 2aa^% k ) 



+ P b a (i + \sy >f t . 



= x 2 



t \"V3 



V N\jj N\i T N\ij N\j T * N\i J ~ 



(3.75) 



Finally, the equation (3.73) becomes the following simple form in terms of the Lagrangian coordinates: 



A' 



t \ -V3 



$ (i) ,^(i) $ (i) + *( 2 )Wo(A 3 



(3.76) 



We are now in the position to be able to solve the shift vector iteratively. The lowest order part of this equation yields 



Ax/3- 



(3)t 



0. 



(3.77) 



so the solution of the shift vector can be chosen as 
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4 3)I = + O (A 2 ) . 

Accordingly, we can conclude that the shift vector /2jf is of the second order and obtain the following equation at the second 
order A 2 of Eq.(^7c|: 



Axffi' = A 2 



(1) 1/3 (-^W q,W + *( 2 ) 

\tlJ V N \n N\i N\ij N\] NH 



(3.78) 



The consistency of our formulation can be easily checked if one confirms that the right-hand side of this equation satisfies 
also divergenceless condition with respect to the Lagrangian coordinates by using Eq.( 3.661). Using Green function of the 
Laplacian, we can express the solution of this equation as 



/3^ ]i (X,t) = A : 
where 



4tt|X -Y\ 



(3.79) 



(3.80) 



This is the expression of the shift vector that we have looked for in order to solve Eq.(3.71). Thus we are able to express all 
perturbation quantities of the metric in our coordinates in terms of the trajectory field. This must be so because it is only 
a dynamical variable in the Lagrangian description. If noting that we have considered only the shift vector (5 r not fli, we 
find that the second order solution of the metric component got has a growing mode: = a 2 f3 l tx t. This time dependence 
exactly agrees with the results derived recently by Matarrese, Mollerach & Bruni ( 199S| ). (The comparison needs the caution 
that they used the conformal coordinates.) However, our approach is essentially different from them: they derived the second 
order solutions of metric perturbations in the Poisson gauge, which reduces to the longitudinal gauge at the lowest order 
with respect to the scalar mode, by transforming the known solutions in the synchronous comoving coordinates according 
to the second order gauge transformation. Similarly, it will be easy to check whether other second order solutions of metric 
perturbations derived by our formalism agree with their results. Here, we are interested in the PN solution of the trajectory 
field. 



By inserting Eq.( 3.79 ) into Eq.( 3.71 ), we can obtain the following equation at the lowest order of A: 

(X,t) = 0. 



dt 2 a dt 

As an irrotational case, the form of this equation allows us to seek the solutions of the following form 
WsP T (X,t) = Q T At)t pN >~f\x), Q T z {ti) = Vx ■ (PA °H (1) = 0, 

and from the above equations we can obtain 

0, 

a 
with 

(FJV) E (1) (x) = n(x), Vx-n = o, 



(3.81) 



(3.82) 



(3.83) 



where II is the unknown function determined by the initial conditions. The solutions of equation ( p.82| ) can be easily found 
to have only a decaying mode. The first order solution may be safely ignored because the decaying mode play physically no 
important role: 



T(PiV) = (l) 



0. 



Next, we similarly consider the second order solution of the PN transverse part with the form 



(PN) S (2)T 



(x,t) = QLW (FJV) sf ) (x) 



(PAT)„(2) 



0. 



(3.84) 



(3.85) 



Inserting this ansatz, Eqs.( 3.7£ ) and ( 3.84 ) into Eq.( 3.71 ), we can obtain the following equation at the second order in A: 



?22 "T" 2 Q Z 2 

a 



Vx x 



(PA0„(2) 



-4/3 



Vx 



The solution of this equation can be found by solving the linear ordinary differential equations: 



(3.86) 



Qzz 2 Qzz 

a 



tj iti) 



-4/3 



(3.87) 
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with 



A (i>JV)g(2) 



(1 



A x f3 T 



(3)i 



{l + Z!)H] 



(1) 



(1) 



N\jj N\i 



(3.88) 



Since we expect the quantity Q zz { to be much smaller than the order A quantities such as 5 at the initial time, we can 

probably adopt again the following initial condition: 



QL(ti) (PJV H (2) 



0. 



Then we can find the solution of Eq.(3.87): 



(1M 



2/3 



12 + 8 



-1/3 



(3.89) 



(3.90) 



Thus the transverse part of the first PN displacement vector has a growing mode even if the Newtonian trajectory field does 
not have a vor ticity flow. We also note that the solution with the growing mode is a particular solution of the differential 
equation ( 3.86 ). From the observational viewpoint, it will be important that we was able to connect the solutions of vector 
perturbations with the initial density fluctuation field which is available through observations of the anisotropies of CMB and 
the two-point correlation function of galaxies. 



4 RESULTS AND IMPLICATIONS 



Summarizing the results derived until the previous section, we have obtained the following family of trajectories as a physical 
complete PN solution up to the second order of A superposed on an isotropic homogeneous deformation solution of Einstein-de 
Sitter background: 



f(X,t) = X + q z (t)V x <f^(X) + q„(t)V x 9$(X) + ±Q%(t)VxV P N(X) 



with 



(!) (h 



t \ 2 / 3 



/3\ 2 f 3 /i\ 4/3 3/t\ 2 / 3 1 4 /t\- 1 
9 " = UJ flit J ~2 + 35 (tlJ 

<*-(§)[§(£)"■ 



2 4 

3 + 15 



J i 



12 + 8 



(r ,s 



and 



Ax* 



9 m ^(i) 

* N\ii* N\jj 
1 



* N\ij w N\ij> 



2 5, 



2 V w l' 



2 15 
+ — 



(*^.(y)*w ( y) 



~ (i + ^) 2 i? 

(1 + Zj )2tf \ w lM ^1' W W lj 



4tt|X 



'3 (2) (X) 



(2) 



(4.1a) 

(4.1b) 
(4.1c) 
(4.1d) 
(4.1e) 

(4.1f) 
(4-lg) 

(4.1h) 
(4.1i) 



We remark that both longitudinal and transverse perturbation solutions of the first PN displacement vector have a growing 
mode and the same time dependence, namely, oc t 2 ^ 3 . In particular, the appearance of the transverse part in PN displacement 
vector is interesting because there is no co rresp ondence in Newtonian order up to the third order in A (Buchertl994 ). The 
physical interpretation is examined by us ( 1998 ) in detail. According to the study, it has been shown that the existence of 
the growing transverse mode guarantees the PN conservation law of vorticity field along flow lines derived by extending the 
Newtonian Kelvin's circulation theorem for the collisionless fluid. 

To see the effect of the PN corrections clearly, let us make a simple order estimation. For the purpose, consider the initial 
density fluctuation 5 g (i) with characteristic length I in the comoving frame. It should be noted that, due to a(to) = 1, we 
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may regard the length I as the physical scale at the present time to. Therefore, if, for convenience, we consider f l (X,t) which 
represents the velocity field of the fluid elements in the comoving frame, it may be estimated roughly at some time t as follows: 



4.* im H.(i fffl )']+Kiffi(^) (^) 2 +^((T T | 7R ) (X W ) 2 + °(A 3 ,c-% (4.2) 



where we used the following estimation: 
V iV(0 ~<M0 1 ' 



*iv(i) ~ (X(o) 



2 



(PJV) = (2) 



(<*s(o) 

(S s (i))V- (4-3) 



Noting ctj = cto/(l + •z/) 3/ ' 2 , q z ~ (1 + Zr)a and so on, the above equation becomes 

f {l) (t) ~l \a8 g{l) (to)+da(S g(l) (t )) 2 ^ + la (J^j (5 g(l) (t )) 2 + ■ ■ ■ , (4.4) 

where S g m(to) is defined by (1 + Zi) S g m representing the present value of the density contrast for scale I extrapolated by the 
linear theory. Note that both the longitudinal and the transverse parts of the PN displacement vector are of the same order 
of magnitude. 

This expression suggests the following interpretations. First, we can conclude that the first and second terms in the 
right-hand side of Eq.( |4.4| ) agree with the result calculated by Newtonian theory if they are expressed in the terms of the 
gauge invariant quantities at the initial time, and the third term is the first nonvanishing PN correction whose order is 
(Z/(cio)) 2 (5 9 (;) (to)) 2 - Note that cto is the horizon scale at present time to- Comparing the first term with the third term on the 
right-hand side, it is also obvious that the expansion parameter of our scheme is the perturbation quantity (l/(cto)) 2 S g ^ (to). 
Consequently, the above expression for the field of trajectories may describe only the congruence of the fluid elements inside 
the scale satisfying the following condition: 

( 5 s(i)(*o)) 1/2 



If applying the usual power spectrum to the above density fluctuation field, this condition means that the solution (4.4) can 
describe the motion of fluid elements up to quasi nonlinear regime for fluctuations not only inside but also beyond the present 
Hubble horizon scale cto, because the magnitude of density contrast at the scale cto is much smaller than unity. The COBE 
microwave background measurement suggests that the power spectrum has the positive slopes on the large scale, supporting 
the assumption of horizon-scale homogeneity on the horizon scale. Furthermore, Eq.( [4.4| ) means that as long as we consider 
the evolution of fluctuation field whose characteristic scale I is much smaller than the horizon scale cto, the fluctuation is well 



described by Newtonian theory because Newtonian correction term, namely the second term in Eq.(4.4), keeps being much 
larger than the first PN correction term from the initial time until today. This may not be correct when one consider the 
fluctuations with larger scale. In fact we can evaluate the redshift z n when the Newtonian correction term begins to become 
larger than the PN correction term for the density fluctuation with the scale I: 

l + z N ~{^-)\ (4.6) 

This redshift is about the same as the redshift when the fluctuation of the scale I just entered into the Hubble horizon ctjv. 
In other words, the PN correction might play an important role for the early evolution of the density contrast. For example, 
SDSS is expected to survey the region of several hundred megaparsecs cube, then the above estimate gives zn ~ 40 for the 
fluctuation with the scale 300/i -1 Mpc. Since the density field at such a redshift is still in the linear stage, the PN effects may 
not be so important for such scales. If this is the case, it will confirm the validity of the use of Newtonian simulation for 
such scales. However, it should be noted that we could come to these conclusions by formulating the Lagrangian perturbation 
theory in the longitudinal gauge and expressing the solution of the trajectory field in terms of the gauge invariant quantities at 
the initial time. These conclusions may be clearer if one recalls that set of equations in the relativistic linearized theory can be 
reduced to a Newtonian-like system of equations expressed in terms of the gauge invariant variables in the longitudinal gauge 
(see Appendix A). Furthermore, this gauge condition allows us to have a straightforward Newtonian limit on scales much 
smaller than the horizon scale even in the nonlinear situations. The PN corrections are small compared with the Newtonian 
solutions. However, because the first PN solutions has a transverse part with a growing mode which cannot be predicted 
by Newtonian theory as well as by gauge invariant linearized theory, its characteristic effects might turn out on the large 
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scales. Actually, if we compared the order of PN transverse part with Newtonian transverse one at third order of A derived 
by Buchert ( |l99^ ), we can estimate the redshift Z(t)jv when the two become the same order of magnitude for the fluctuation 



with the scale 50/i 1 Mpc as 

l + z {T)N ~4. (4.7) 
This shows that the PN transverse mode might play an important role in the observational cosmology. In particular, the effects 



may induce the secondary temperature anisotropics of CMB as the Rees-Sciama effect (1968) and potentially be observable 



Mollerach & Matarrese ( 1997 ) investigated the secondary CMB anisotropics induced by the second order metric perturbations 



during the travel of CMB photons from the last scattering surface to us, using the formalism developed by Pyne and Carroll 



(1996). They concluded that the effect of the higher order perturbations is much smaller than that of the usual gravitational 
lensing due to the first order scalar perturbation. However, the study of second order anisotropies will be relevant because 
they have a possibility to produce a non negligible contribution compared to the first order ones, due to the long distances 
of the travel. The reason is that, since the second order effects are caused by the integrals of the non-linearities in the metric 
perturbations along the photon path, they accumulate the small effects. Moreover, the second order effects do not produce the 
cancellation for the integral along photons path in contrast with the linear effects, so they may give the primary contribution 
to some statistical measures as, for example, the three-point function of temperature anisotropies. Thus, it is important to 
estimate the magnitude of anisotropies due to second order vector perturbations as a possible contribution to the linear order 
anisotropy calculations. It should be emphasized that our formalism allows us to estimate the magnitude of the second order 
vector mode of metric perturbations, namely the shift vector, by using the power spectrum of density fluctuation field. We 



are now investigating the secondary effect using the "total angular momentum method' developed by Hu & White (1997). 
According to the formalism, it seems that the growing vector perturbation produces a behavior of anisotropy largely different 
form the scalar mode results. It will be interesting to compare our result with that of Mollerach & Matarrese. Furthermore, for 
the N-body simulations of more wide region of the universe, the PN corrections might also play important role. Our approach 
has an advantage to be used in numerical work based on the presently available Newtonian simulation. In any case, it would 
be very interesting to see these effects due to the PN corrections. These works are now in progress. 



5 SUMMARY 

In this paper, we have formulated the PN Lagrangian perturbation theory in the perturbed Einstein-de Sitter universe using 
the (3+1) formalism and investigated their effects on the evolution of the large-scale structure of the universe up to the quasi 
nonlinear regime. In this formalism, all force terms in the first PN equations of motion are expressed in terms of the (gauge 
invariant) Newtonian quantities. As a result, we can solve the trajectory field iteratively in the PN approximation. 



Our formulation based on the gauge condition (2.26) has a natural Newtonian limit and, consequently, we can show that 



(3) 

the longitudinal part of shift vector p L and slice condition 8, which are freely specified within the residual gauge freedom, 
do not influence the solution of the trajectory field only if the Newtonian one is assumed to be irrotational. 

It has been shown that the longitudinal solutions of the first PN displacement vector in our coordinate conditions have 
spurious modes caused by the gauge ambiguities and, if the solution is expressed in the terms of the gauge invariant quantities 
at the initial time, the PN irrotational lowest order solution (order A) take exactly same expression with corresponding 
Newtonian first order solution. This allows us to regard the gauge invariant lowest order solution as a physical Newtonian 
lowest order one. On the other hand, since the transverse solution of the PN displacement vector is generated from the 
beginning by the gauge invariant quantity, namely, the transverse part of the shift vector /3y , it first appears with the same 
order of the first non- vanishing PN irrotational solution (order A 2 ) expressed in the terms of the gauge invariant quantities at 
the initial time. Then it was shown that the perturbation quantity (l/(cto)) 2 S g (to) become an essential expansion parameter in 
the PN approximation in the case that we consider the evolution of the fluctuations in the expanding universe. The smallness 
of parameter will explain why Newtonian treatment in our coordinate choice is so good approximation even for horizon scale 
fluctuations. This consideration seems to be always the case if we study the behavior of nonrelativistic matter in the coordinate 
which has a straightforward Newtonian limit. 

We have explicitly derived physical nonvanishing PN corrections to the Newtonian dynamics and found some interesting 
results. One is that the first PN displacement vector has a transverse part with a growing mode, which is not present in 
cosmological Newtonian theory as well as in the gauge invariant linearized theory. The effect may be potentially observed as 
a characteristic pattern of the large-scale structure. The other is that we was able to estimate rigorously the PN corrections 
to the Newtonian approximation in the situation we are interested in. 

Since the observable region of the large-scale structure increases steadily, we expect that our work will play some important 
roles in the future. 
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APPENDIX A: GAUGE INVARIANT LINEARIZED THEORY 

In this appendix, we consider the linearized theory in the Einstein-de Sitter universe (i.e., K. — A = 0) to determine the 
behavior of the metric quantities and the matter variables in the early stage of the evolution of the density fluctuation. In 
the linearized theory, all perturbation quantities can be divided into the scalar parts, the vector parts and tensor parts with 
respect to the background geometry: 

(3' =p L<i + fa, 

V 1 = VL,i + Vt, 

(Al) 

where 

0T,i = Vr,i — 0, ■ ■ ■ , 



and it should be noted that we consider only /3 ! not /3j. The equations for scalar modes under our gauge conditions (2.26) 
become 

1 dv L 2 a 1 d(i L 2 d a 1 . , . „, 

c at cr a c at c a cr 

3 a _ 3 dip . „ , . 

0= — a+ + A f p L , A4 
c a c at 

- 5 V-^( { + 3a)-- I ---e, (A5) 

^ Afqp= ^ d + i±e, (A6) 
a 2 c 2 c a 



dip & _ 2 /l 1 

at a 



AirGp b a (-5VL + -Pt^ , (A7) 
ire a = a — 1. From the above equations, the equation for S can be derived as 

9 +2^-4»rGM>=0. (A8) 



at 2 a at 

where 

o s = — jaa 



Here, <5 9 denotes a gauge invariant quantity, and from Eq.(A£), we obtain the two evolution modes of the density fluctuation 



S g oc a and oc a 3 / 2 . It should be remarked that this quantity is the only physical density fluctuation field. We can make use 
of the residual gauge freedom for Pl'- 

Pl = 0. (A10) 
Then the above equations become 

A/a = ^-8 g , All 

c 2 

AnGa 2 p b ,.,„, 

A /^ = ^ 5 s> ( A12 ) 



The Eqs.(^il[) and (|m|) entirely agree with the Newtonian-like Poisson's equation. However, since general relativity doesn' 



t 

depend on the scale, Eq.(All) is valid for the density fluctuation of all scale in the linear regime. It is noted that the structure 
of set of the above equations is just like Newtonian theory but not exactly same because the gauge invariant linearized theory 
must be constrained by the condition i 9 C 1. Thus we can conclude that the solutions for the expansion 

a = —5- H j- H , 

c 2 c 4 
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c* c 4 



become 



A/Q (2) = ^Ga Pb 8 g , 



„( 4 ) 



,(«) 



0, 



(A14) 
(A15) 

A/V^ 2 ' = ATrGa 2 p b 5 g , (A16) 
V> (4) = V> (6) = • • • = 0. (A17) 

Thus if we consider only the growing mode of 8 g , we immediately find that the solution for becomes oc a . Then the S 
under our gauge conditions becomes from Eqs.(A13) and (A.9) 

v( 2 ) 



5 + 2- 



(A18) 



The 5 has the spurious mode 2q (2) /c 2 . This means that if we use 5 instead of 5 g , the undesirable mode will appear. Therefore, 
in order to see the physical density fluctuation, we must choose 5 g as a density fluctuation instead of 5. 

Finally, we give a physical meaning of the gauge condition (A10). It is easily seen that under the gauge condition the 
extrinsic curvature which represents the cosmic isotropic expansion becomes from Eqs.(2.15) and (A4) 

1 da(r) 1 da 



K 



a dr 



a(l + a) cdt' 



(A19) 



where we used the ip remains constant from Eq.(A12) in the linear regime. Thus the gauge condition corresponds to the choice 
of the proper time r on the spatial hypersurface. 
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